home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Turnbull China Bikeride
/
Turnbull China Bikeride - Disc 1.iso
/
ARGONET
/
PD
/
MATHS
/
RLAB
/
RLAB125.ZIP
/
!RLaB
/
toolbox
/
fmin
< prev
next >
Wrap
Text File
|
1994-09-20
|
3KB
|
124 lines
//---------------------------------------------------------------------------
// fmin.r:
// x = fmin(f,ax,bx,tol,print)
//
// find the minimum of the function f(x) in the interval [ax,bx]
//
// f = The real-values function of a single variable.
//
// ax,bx = left and right endpoints of search interval
//
// tol = maximum length of interval of uncertainty of minimum
// tol >= 0
//
// print = an optional last agrument, if nonzero a printed tracing of
// the steps is given
//
// x = the value of x approximating where f attains a minimum in [ax,bx]
// This algorithm is from: "Computer Methods for Mathematical
// Computations", Forsythe, Malcolm, and Moler, Prentice-Hall, 1976.
// Adapted by: Duane Hanselman, University of Maine, (207)-581-2246
// initialization
//---------------------------------------------------------------------------
fmin = function ( f, ax, bx, tol, _print)
{
global (eps)
if (!exist (_print)) { print = 0; else print = _print; }
num = 1;
seps = sqrt(eps);
c = 0.5*(3.0 - sqrt(5.0));
a = ax; b = bx;
v = a + c*(b-a);
w = v; x = v;
d = 0.0; e = 0.0;
fx = f (x);
if (print) { fmin_data = [1, x, fx] ? }
fv = fx; fw = fx;
xm = 0.5*(a+b);
tol1 = seps*abs(x) + tol/3.0;
tol2 = 2.0*tol1;
// main loop
while (abs (x-xm) > (tol2 - 0.5*(b-a)))
{
num = num+1;
gs = 1;
// is parabolic fit possible
if (abs(e) > tol1) // yes, so fit parabola
{
gs = 0;
r = (x-w)*(fx-fv);
q = (x-v)*(fx-fw);
p = (x-v)*q-(x-w)*r;
q = 2.0*(q-r);
if (q > 0.0) { p = -p; }
q = abs (q);
r = e;
e = d;
// is parabola acceptable
if ((abs (p) < abs (0.5*q*r)) && (p > q*(a-x)) && (p<q*(b-x)))
{
// yes, parabolic interpolation step
d = p/q;
u = x+d;
step = " num x fx parabolic";
// f must not be evaluated too close to ax or bx
if (((u-a) < tol2) || ((b-u) < tol2))
{
si = sign(xm-x) + ((xm-x) == 0);
d = tol1*si;
}
else // not acceptable, must do a golden section step
gs=1;
}
}
if (gs) // a golden-section step is required
{
if (x >= xm) { e = a-x; else e = b-x; }
d = c*e;
step = " num x fx golden ";
}
// f must not be evaluated too close to x
si = sign(d) + (d == 0);
u = x + si * max ([abs (d), tol1]);
fu = f (u);
if (print) { fmin_data = [num, u, fu] ? disp(step) ? }
// update a, b, v, w, x, xm, tol1, tol2
if (fu <= fx)
{
if (u >= x) { a = x; else b = x; }
v = w; fv = fw;
w = x; fw = fx;
x = u; fx = fu;
else // fu > fx
if (u < x) { a = u; else b = u; }
if ( (fu <= fw) || (w == x) )
{
v = w; fv = fw;
w = u; fw = fu;
else if ( (fu <= fv) || (v == x) || (v == w) ) {
v = u; fv = fu;
}}
}
xm = 0.5*(a+b);
tol1 = seps*abs(x) + tol/3.0;
tol2 = 2.0*tol1;
}
return x
};